library(readxl)
data <- read_excel("Dataset3.xlsx")
head(data)
## # A tibble: 6 × 4
## Viti ProdhimiDritherave PlehraKimike CO2
## <dbl> <dbl> <dbl> <dbl>
## 1 1961 293932 6650 1.29
## 2 1962 313508 4110 1.35
## 3 1963 293449 4040 1.10
## 4 1964 350144 3990 1.03
## 5 1965 335293 3830 1.08
## 6 1966 433094 5560 1.23
Ky projekt synon të analizojë ndryshimet dhe trendet e prodhimit të drithërave në lidhje me përdorimin e plehrave kimikë dhe emëtimet e CO2 në një periudhë afatgjatë, duke përdorur analiza të serive kohore.
library(TTR)
# Moving Average (MA) me periudhë 3
ma_data <- ma(ts_prodhimi_dritherave, order = 3)
plot(ts_prodhimi_dritherave, main = "Seri Kohore dhe Moving Average", col = "blue", type = "l")
lines(ma_data, col = "green", lwd = 2)
legend("topright", legend = c("Seri Kohore", "Moving Average"), col = c("blue", "green"), lwd = 2)
Ky grafik paraqet serinë kohore të prodhimit të drithërave (1961-2020) dhe mesataren e lëvizshme të saj me një periudhë prej 3 vitesh.
Vijat blu përfaqësojnë serinë kohore origjinale, ndërsa vija jeshile paraqet mesataren e lëvizshme.
Vijat blu tregojnë luhatjet e forta në prodhimin e drithërave, duke përfshirë disa luhatje të papritura, sidomos për vitet e hershme.
Mesatarja e lëvizshme (vija jeshile) zbut këto luhatje, duke dhënë një trend të qartë rritës të përgjithshëm gjatë periudhës.
Përdorimi i mesatares së lëvizshme është i dobishëm për të identifikuar trendin afatgjatë duke eliminuar luhatjet e rastësishme.
Shohim që vijat përshtaten mirë me njëra-tjetrën.
4.2 Mesatarja e lëvizshme eksponenciale për prodhimin e drithërave:
ema_data <- EMA(ts_prodhimi_dritherave, n = 3)
plot(ts_prodhimi_dritherave, main = "Seri Kohore dhe Exponential Moving Average", col = "blue", type = "l")
lines(ema_data, col = "purple", lwd = 2)
legend("topright", legend = c("Seri Kohore", "Exponential Moving Average"), col = c("blue", "purple"), lwd = 2)
EMA ka një peshë më të madhe për vlerat më të fundit të serisë kohore, duke reaguar më shpejt ndaj ndryshimeve sesa mesatarja e lëvizshme.
Vijat blu tregojnë serinë origjinale, ndërsa vijat vjollcë përfaqësojnë EMA.
EMA përshtatet më mirë me ndryshimet e fundit të prodhimit të drithërave. Për shembull, reagimi ndaj rritjes së shpejtë pas viteve 1990 është më i shpejtë në vijën vjollcë.
Është më e ndjeshme ndaj ndryshimeve të kohëve të fundit në krahasim me MA.
EMA është një mjet më i mirë kur është e rëndësishme të kapet sjellja e fundit e serisë kohore.
4.3 Sheshimi i thjeshtë eksponencial për prodhimin e drithërave:
data_ts <- ts(data$ProdhimiDritherave, start = c(1961), frequency = 1)
shesh_exp_drithera <- HoltWinters(ts_prodhimi_dritherave, beta = FALSE, gamma = FALSE)
shesh_exp_drithera
## Holt-Winters exponential smoothing without trend and without seasonal component.
##
## Call:
## HoltWinters(x = ts_prodhimi_dritherave, beta = FALSE, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 0.9999202
## beta : FALSE
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 684021.6
plot(shesh_exp_drithera,
main = "Sheshimi i thjeshtë eksponencial për prodhimin e drithërave",
xlab = "Viti",
ylab = "Prodhimi i drithërave")
Ky grafik paraqet një model SES pa trend dhe sezonalitet.
Linjat e parashikimit tregojnë një mesatare të zbutur të serisë kohore.
Modeli SES përshtatet mirë me të dhënat, duke kapur trendin rritës në seri. Vlerat e alfa-s sugjerojnë një peshim të fortë të vlerave të fundit.
Nuk mund të trajtojë trendet ose sezonalitetet afatgjata, por është efektiv për seritë pa këto karakteristika.
SES është një metodë e thjeshtë për zbutjen e të dhënave dhe përshtatet mirë për të dhëna me variacion minimal.
4.4 Sheshimi eksponencial me Holt-Winters për prodhimin e drithërave:
# Aplikimi i modelit Holt-Winters me alpha dhe beta (pa gamma)
holt_model <- HoltWinters(data_ts, beta = 0.1, gamma = FALSE) # Gamma = FALSE heq sezonalitetin
# Shfaqja e parametrave të modelit
holt_model
## Holt-Winters exponential smoothing with trend and without seasonal component.
##
## Call:
## HoltWinters(x = data_ts, beta = 0.1, gamma = FALSE)
##
## Smoothing parameters:
## alpha: 0.9999525
## beta : 0.1
## gamma: FALSE
##
## Coefficients:
## [,1]
## a 684022.202
## b 2815.481
# Vizualizimi i modelit Holt-Winters
plot(holt_model,
main = "Modeli Holt-Winters për Prodhimin e Drithërave (Me Trend)",
xlab = "Viti",
ylab = "Prodhimi i Drithërave")
Ky model përfshin trendin në mënyrë të qartë, por nuk përfshin sezonalitet. - Vijat e modeluara tregojnë një trend të lëmuar të rritjes së prodhimit të drithërave. - Modeli është shumë i përshtatshëm për të dhënat që tregojnë një trend rritës afatgjatë, si në këtë rast. - Parametrat alfa dhe beta përputhen mirë me të dhënat origjinale, duke rezultuar në një model të saktë. - Konkluzion: - Holt-Winters është një metodë më e avancuar për të trajtuar seri me trend dhe siguron një përshtatje më të mirë krahasuar me SES.
library(forecast)
# Modeli ARIMA për Prodhimi i Drithërave
fit_eme_CO2 <- auto.arima(ts_CO2)
summary(fit_eme_CO2)
## Series: ts_CO2
## ARIMA(1,0,0) with non-zero mean
##
## Coefficients:
## ar1 mean
## 0.8695 1.5789
## s.e. 0.0583 0.2436
##
## sigma^2 = 0.07664: log likelihood = -7.76
## AIC=21.53 AICc=21.96 BIC=27.81
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE
## Training set 0.006669262 0.2721802 0.1926131 -3.439636 13.6033 1.01612
## ACF1
## Training set -0.05593145
# Parashikimi me ARIMA
forecast_eme_CO2 <- forecast(fit_eme_CO2, h = 4)
autoplot(forecast_eme_CO2)
Ky grafik tregon modelin ARIMA(1,0,0) për emetimet e CO2, duke përfshirë vijën e parashikimit dhe kufijtë e besueshmërisë. - Modeli përputhet mirë me të dhënat origjinale të emetimeve të CO2. Autokorrelacionet e vonesave (lags) tregojnë se një urdhër i parë autoregresiv (AR) është i mjaftueshëm. - Kufijtë e besueshmërisë janë të ngushtë, duke sugjeruar një parashikim të besueshëm. - ARIMA është një metodë e fuqishme për të trajtuar seri kohore komplekse. Për CO2, modeli me komponent autoregresiv është i mjaftueshëm për kapjen e dinamikës së të dhënave.
library(forecast)
# Modeli ARIMA për Prodhimi i Drithërave
fit_plehra <- auto.arima(ts_plhera_kimike)
summary(fit_plehra)
## Series: ts_plhera_kimike
## ARIMA(1,0,0) with non-zero mean
##
## Coefficients:
## ar1 mean
## 0.9063 44410.14
## s.e. 0.0505 16696.52
##
## sigma^2 = 196055624: log likelihood = -657.8
## AIC=1321.6 AICc=1322.02 BIC=1327.88
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 938.7646 13766.64 8955.413 -19.40321 37.64379 1.07692 0.06175204
# Parashikimi me ARIMA
forecast_plehra <- forecast(fit_plehra, h = 4)
autoplot(forecast_plehra)
Ky grafik paraqet modelin ARIMA(1,0,0) për sasinë e plehrave kimike, së bashku me parashikimet për vitet e ardhshme. - Modeli kap rritjen e qëndrueshme të plehrave kimike gjatë periudhës. Autokorrelacionet sugjerojnë një trend të fortë autoregresiv për këtë variabël. - Vijat e parashikimeve tregojnë një trend rritës të vazhdueshëm, që është në përputhje me trendin historik. - ARIMA(1,0,0) është efektiv për trajtimin e kësaj serie, duke ofruar parashikime të besueshme për përdorimin e plehrave kimike.
data_train = head(data, 48)
data_test = tail(data, 12)
5.2.1) Ndërtimi i modelit dhe saktësia e tij
library(forecast)
library(fpp2)
## ── Attaching packages ────────────────────────────────────────────── fpp2 2.5 ──
## ✔ fma 2.5 ✔ expsmooth 2.3
##
drift.drithera <- rwf(data_train$ProdhimiDritherave, h = 12, drift = TRUE)
accuracy(drift.drithera, data_test$ProdhimiDritherave)
## ME RMSE MAE MPE MAPE MASE
## Training set -1.888653e-11 95967.30 57358.58 -0.9262221 9.728369 0.9755229
## Test set 3.549491e+04 45834.02 38969.84 5.0902108 5.608772 0.6627773
## ACF1
## Training set 0.1057804
## Test set NA
5.2.2) Kontrolli i mbetjeve të modelit
checkresiduals(drift.drithera)
##
## Ljung-Box test
##
## data: Residuals from Random walk with drift
## Q* = 8.081, df = 10, p-value = 0.6209
##
## Model df: 0. Total lags used: 10
Rezultatet nga grafiku i mbetjeve: - Grafiku i mbetjeve tregon se ato janë shpërndarë rastësisht rreth vijës zero, që është një tregues pozitiv se modeli nuk ka mbetje të strukturuara.
Autokorrelacioni i mbetjeve: - ACF për mbetjet tregon mungesë të
autokorrelacioneve domethënëse, që sugjeron se mbetjet janë në mënyrë të
pavarur nga njëra-tjetra dhe janë stacionare.
Testi Ljung-Box: - Rezultati i testit Ljung-Box tregon një p-vlerë prej
0.6209, që është mbi pragun 0.05. Kjo sugjeron se nuk ka autokorrelacion
të rëndësishëm mes mbetjeve.
5.3.1) Krijimi i modelit dhe rezultatet e tij
# Ngarkimi i bibliotekave
library(forecast)
library(ggpubr)
##
## Attaching package: 'ggpubr'
## The following object is masked from 'package:forecast':
##
## gghistogram
# Krijimi i serisë kohore për trajnimin
ts_prodhimi_dritherave <- ts(data_train$ProdhimiDritherave, start = 1961, frequency = 1)
# Krijimi i modelit ARIMA me auto.arima
arima.drithera <- auto.arima(ts_prodhimi_dritherave)
# Rezultatet e modelit
summary(arima.drithera)
## Series: ts_prodhimi_dritherave
## ARIMA(0,1,0)
##
## sigma^2 = 9.255e+09: log likelihood = -605.98
## AIC=1213.95 AICc=1214.04 BIC=1215.8
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 6559.624 95193.06 57578.95 0.2119825 9.723723 0.9792708 0.1055468
5.3.2)Parashikimi me modelin ARIMA dhe paraqitja grafike
# Parashikimi me modelin ARIMA
forecast_arima <- forecast(arima.drithera, h = 12)
# Vizualizimi i parashikimit
plot(forecast_arima, main = "Parashikimi me Modelin ARIMA", ylab = "Prodhimi i Drithërave", xlab = "Viti")
lines(ts(data_test$ProdhimiDritherave, start = 2009, frequency = 1), col = "red", lwd = 2)
lines(ts(arima.drithera$fitted, start = 1961, frequency = 1), col = "orange", lwd = 2)
legend("topleft", c("Train", "Fit", "Prediction", "Test"),
fill = c("black", "orange", "blue", "red"), box.col = "black")
- Kufijtë e besueshmërisë janë relativisht të ngushtë, duke treguar
besueshmëri të lartë në parashikime. - Të dhënat reale dhe parashikimet
ndjekin një trend të ngjashëm, megjithëse ka disa devijime në vlerat
specifike.
5.3.3) Saktësia e modelit ARIMA
accuracy(forecast_arima, data_test$ProdhimiDritherave)
## ME RMSE MAE MPE MAPE MASE
## Training set 6559.624 95193.06 57578.95 0.2119825 9.723723 0.9792708
## Test set 78999.000 81525.74 78999.00 11.4104076 11.410408 1.3435711
## ACF1
## Training set 0.1055468
## Test set NA
Performanca e modelit është më e mirë për setin e trajnimit se për setin e testimit, me disa metrika si RMSE, MAE, dhe MAPE që janë më të larta për setin e testimit, duke sugjeruar se modeli nuk është shumë i saktë në parashikimet për të dhënat që nuk ka parë më parë (setin e testimit).
5.3.4) Kontrolli i mbetjeve të modelit ARIMA
checkresiduals(arima.drithera)
##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,0)
## Q* = 8.1954, df = 10, p-value = 0.6098
##
## Model df: 0. Total lags used: 10
# Testi i Ljung-Box për autocorrelacionin e mbetjeve
Box.test(arima.drithera$residuals, lag = 10, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: arima.drithera$residuals
## X-squared = 8.1954, df = 10, p-value = 0.6098
Veme re se autokorrelacionet e mbetjeve jane brenda intervalit te besimit dhe kjo do te thote se seria e mbetjeve eshte stacionare. Rezultatet sugjerojnë që nuk ka autokorrelacion të rëndësishme në mbetjet e modelit. P-vlera prej 0.6098 është më e madhe se 0.05, që nënkupton se mbetjet janë të pavarura dhe modelin ARIMA(0,1,0) mund të jetë i përshtatshëm për të dhënat.
5.3.5) Vizualizimi i shpërndarjes së mbetjeve
ggdensity(as.numeric(arima.drithera$residuals), fill = "lightgray", main = "Shpërndarja e Mbetjeve")
ggqqplot(as.numeric(arima.drithera$residuals), color = "magenta",main = "QQ-Plot i Mbetjeve")
Që modeli te jetë i mirë ne presim që kuantilet teorike dhe të zgjedhjes të jenë brenda zonës gri pra afer drejtëzës. Në këtë rast shumica e tyre janë brenda zonës gri.
5.4.1) Krijimi i modelit RegARIMA dhe rezultatet e tij
# Ngarkimi i bibliotekave të nevojshme
library(forecast)
library(ggpubr)
library(tseries)
# Krijimi i regresorit si seri kohore (PlehraKimike)
ts_plehra <- ts(data_train$PlehraKimike, start = 1961, frequency = 1)
# Krijimi i modelit RegARIMA me regresorë PlehraKimike
regarima_model <- auto.arima(
ts_prodhimi_dritherave,
xreg = cbind(ts_plehra),
seasonal = FALSE
)
# Rezultatet e modelit
summary(regarima_model)
## Series: ts_prodhimi_dritherave
## Regression with ARIMA(1,0,0) errors
##
## Coefficients:
## ar1 intercept xreg
## 0.7611 408293.30 4.7118
## s.e. 0.1053 50996.14 0.7259
##
## sigma^2 = 5.179e+09: log likelihood = -603.82
## AIC=1215.64 AICc=1216.57 BIC=1223.13
##
## Training set error measures:
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 3379.92 69679.09 53657.94 -0.937631 9.133682 0.9125844 0.06813269
Ky është një model regresioni që përfshin gabime të ARIMA(1,0,0), dhe ka pasur performancë të mirë për setin e trajnimit, me metrika të sakta si ME, RMSE, dhe MAE që janë relativisht të vogla. AIC, AICc dhe BIC sugjerojnë se ky model mund të jetë një zgjedhje e mirë për këtë dataset.
5.4.2) Parashikimi dhe vizualizimi i tyre
# Regresorët për testimin
xreg_test <- cbind(data_test$PlehraKimike)
# Parashikimi për 12 vitet e ardhshme
forecast_regarima <- forecast(regarima_model, xreg = xreg_test, h = 12)
# Vizualizimi i parashikimeve
plot(forecast_regarima, main = "Parashikimi me Modelin RegARIMA",
ylab = "Prodhimi i Drithërave", xlab = "Viti")
lines(ts(data_test$ProdhimiDritherave, start = 2009, frequency = 1), col = "red", lwd = 2)
lines(ts(regarima_model$fitted, start = 1961, frequency = 1), col = "orange", lwd = 2)
legend("topleft", c("Train", "Fit", "Prediction", "Test"),
fill = c("black", "orange", "blue", "red"), box.col = "black")
Modeli parashikon mirë vlerat por me disa luhatje të vogla.
5.4.3) Saktësia e modelit RegARIMA
accuracy(forecast_regarima, data_test$ProdhimiDritherave)
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 3379.92 69679.09 53657.94 -0.937631 9.133682 0.9125844 0.06813269
## Test set 44722.91 56753.30 46993.34 6.432737 6.789821 0.7992366 NA
5.4.4) Kontrolli i mbetjeve të modelit
# Kontrollimi i rezidualeve
checkresiduals(regarima_model)
##
## Ljung-Box test
##
## data: Residuals from Regression with ARIMA(1,0,0) errors
## Q* = 9.2685, df = 9, p-value = 0.4129
##
## Model df: 1. Total lags used: 10
# Testi Ljung-Box për autocorrelacionin e mbetjeve
Box.test(regarima_model$residuals, lag = 10, type = "Ljung-Box")
##
## Box-Ljung test
##
## data: regarima_model$residuals
## X-squared = 9.2685, df = 10, p-value = 0.5068
Veme re se autokorrelacionet e mbetjeve jane brenda intervalit te besimit dhe kjo do te thote se seria e mbetjeve eshte stacionare. Rezultatet sugjerojnë që nuk ka autokorrelacion të rëndësishme në mbetjet e modelit. P-vlera prej 0.4129 është më e madhe se 0.05, që nënkupton se mbetjet janë të pavarura dhe modeli mund të jetë i përshtatshëm për të dhënat.
5.4.5) Shpërndarja e mbetjeve
# Shpërndarja e rezidualeve
ggdensity(as.numeric(regarima_model$residuals), fill = "lightblue", main = "Shpërndarja e Rezidualeve")
# QQ-Plot i rezidualeve
ggqqplot(as.numeric(regarima_model$residuals),color = "red" ,main = "QQ-Plot i Rezidualeve")
Që modeli të jetë i mirë ne presim që kuantilet teorike dhe të zgjedhjes të jenë brenda zonës gri pra afer drejtëzës. Në këtë rast shumica e tyre janë brenda zonës gri.
5.4.6) Paraqitja e gabimeve e Regresit dhe ARIMA
cbind("Gabimet e regresit" = residuals(regarima_model, type="regression"),
"Gabimet ARIMA" = residuals(regarima_model, type="innovation")) %>%
autoplot(facets=TRUE)
5.5.1) Ndërtimi i modelit
library(forecastHybrid)
## Loading required package: thief
hyb.mod.8 <- hybridModel(ts(data_train$ProdhimiDritherave,start=1961,frequency=1), models = "an",
a.args = list(xreg = ts(data_train$PlehraKimike,start=1961,frequency=1)),
n.args = list(xreg = ts(data_train$PlehraKimike,start=1961,frequency=1)))
## Fitting the auto.arima model
## Fitting the nnetar model
hyb.f8<-forecast(hyb.mod.8,xreg =ts(data_test$PlehraKimike,start=c(2009,1),frequency=1))
## Warning in forecast.hybridModel(hyb.mod.8, xreg = ts(data_test$PlehraKimike, :
## The number of rows in xreg should match h. Setting h to nrow(xreg).
autoplot(hyb.f8)+ylab("Prodhim Drithera")+
autolayer(ts(data_test$ProdhimiDritherave,start=c(2009,1),frequency = 1),series="Data")
Modeli hibrid ndjek mirë trendin e prodhimit të drithërave gjatë periudhës së trajnimit. Vijat e parashikimit janë të përshtatshme dhe të njëtrajtshme me të dhënat reale. Parashikimet për periudhën e testimit janë të afërta me vlerat reale, duke treguar aftësi të mirë të modelit për të kapur dinamikën e të dhënave.
5.5.2) Rezultatet dhe saktësia e modelit
summary(hyb.mod.8)
## Length Class Mode
## auto.arima 19 forecast_ARIMA list
## nnetar 17 nnetar list
## weights 2 -none- numeric
## frequency 1 -none- numeric
## x 48 ts numeric
## xreg 2 -none- list
## models 2 -none- character
## fitted 48 -none- numeric
## residuals 48 ts numeric
accuracy(hyb.f8)
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 2763.768 57771.71 44621.36 -0.668358 7.425179 0.7588952 0.09308877
Vlerat e ulëta tregojnë një përshtatje të saktë të modelit me të dhënat historike. MAPE (7.38%): Gabimi mesatar relativisht i ulët tregon një model të saktë dhe të besueshëm për parashikimet. Përfshirja e plehrave kimike: Regresori i jashtëm përmirëson parashikimet duke ofruar informacion shtesë mbi ndikimin e jashtëm në prodhimin e drithërave
5.5.3) Kontrolli për mbetjet e modelit
checkresiduals(hyb.mod.8$residuals)
##
## Ljung-Box test
##
## data: Residuals
## Q* = 8.477, df = 10, p-value = 0.5824
##
## Model df: 0. Total lags used: 10
Shpërndarja e mbetjeve: - Mbetjet janë shpërndarë rastësisht rreth vijës zero, duke treguar se modeli nuk ka mbetje të strukturuara.
Autokorrelacioni i mbetjeve: - Grafiku ACF për mbetjet nuk tregon autokorrelacione domethënëse.
Testi Ljung-Box: - P-vlera e testit është 0.5903, që është mbi pragun 0.05, duke treguar mungesë të autokorrelacionit të mbetjeve.
5.5.4) Box Test për mbetjet e modelit
Box.test(x = hyb.mod.8$residuals)
##
## Box-Pierce test
##
## data: hyb.mod.8$residuals
## X-squared = 0.40728, df = 1, p-value = 0.5234
p-value = 0.5129 > 0.05, që do të thotë se nuk ka evidence për autokorrelacion në mbetje. Kjo sugjeron se modeli ka kapur informacionin e nevojshëm nga të dhënat dhe mbetjet janë të pavarura.
5.5.5) Paraqitja grafike e mbetjeve të modelit
library(ggpubr)
# Density plot
ggdensity(hyb.mod.8$residuals, fill = "lightgray")
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_density()`).
# QQplot
ggqqplot(hyb.mod.8$residuals)
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Warning: Removed 1 row containing non-finite outside the scale range
## (`stat_qq()`).
## Warning: Removed 1 row containing non-finite outside the scale range (`stat_qq_line()`).
## Removed 1 row containing non-finite outside the scale range (`stat_qq_line()`).
Mbetjet nuk ndjekin një shpërndarje normale. Që modeli te jetë i mirë ne presim që kuantilet teorike dhe të zgjedhjes të jenë brenda zonës gri pra afer drejtëzës. Në këtë rast shumica e tyre janë brenda zonës gri.
5.5.6) Shapiro Test për mbetjet e modelit
shapiro.test(hyb.mod.8$residuals)
##
## Shapiro-Wilk normality test
##
## data: hyb.mod.8$residuals
## W = 0.95, p-value = 0.04332
Mbetjet nuk kanë shpërndarje normale.
5.6.1) Sheshimi i thjeshtë eksponencial dhe saktësia e tij
model1 <- ses(data_train$ProdhimiDritherave, h = 12)
accuracy(model1, data_test$ProdhimiDritherave)
## ME RMSE MAE MPE MAPE MASE
## Training set 6553.061 95194.08 57574.23 0.2095708 9.722069 0.9791905
## Test set 79010.497 81536.88 79010.50 11.4120814 11.412081 1.3437667
## ACF1
## Training set 0.1056043
## Test set NA
Modeli SES ka një MAPE prej 9.72% për setin e trajnimit dhe një 11.41% për setin e testit, duke treguar një përfundim të ngjashëm për të dyja periudhat. Ky është një model i thjeshtë, por mund të ketë vështirësi për të kapur dinamikën më komplekse të serisë kohore.
5.6.2) Sheshimi eksponencial me 2 parametra (metoda Holt) dhe saktësia e tij
# Sheshimi eksponencial me 2 parametra (metoda Holt)
model2 <- holt(data_train$ProdhimiDritherave, h = 12)
accuracy(model2, data_test$ProdhimiDritherave)
## ME RMSE MAE MPE MAPE MASE
## Training set -13579.06 96733.31 55634.80 -3.006817 9.561125 0.9462058
## Test set 51189.36 56402.54 51189.36 7.370276 7.370276 0.8706002
## ACF1
## Training set 0.07601034
## Test set NA
Modeli Holt ka një përformancë të mirë në setin e testimit (MAPE 7.37%) dhe ka përmirësuar performancën krahasuar me SES. Kjo tregon se modeli është në gjendje të kapë më mirë trendet në të dhëna.
5.6.3) Sheshimi eksponencial (me trend të zbutur) dhe saktësia e tij
# Sheshimi eksponencial (me trend të zbutur)
model3 <- holt(data_train$ProdhimiDritherave, damped = TRUE, h = 12)
accuracy(model3, data_test$ProdhimiDritherave)
## ME RMSE MAE MPE MAPE MASE
## Training set -3198.034 93969.01 54998.30 -1.681557 9.496265 0.9353805
## Test set 73535.286 76158.92 73535.29 10.617319 10.617319 1.2506473
## ACF1
## Training set 0.08086838
## Test set NA
Ky model ka një MAPE të lartë në setin e testimit (10.62%) dhe një RMSE të lartë, që sugjeron se modeli me trend të zbutur mund të jetë i dobishëm për disa lloje të serive kohore, por mund të kërkojë përshtatje të mëtejshme për të kapur më mirë tendencat.
5.6.4) Modeli ETS dhe saktësia e tij
# Modeli ETS
model4 <- forecast(ets(data_train$ProdhimiDritherave), h = 12)
accuracy(model4, data_test$ProdhimiDritherave)
## ME RMSE MAE MPE MAPE MASE
## Training set 8326.474 100393.8 61572.48 0.1655844 10.53745 1.047191
## Test set 110880.276 112694.5 110880.28 16.0519018 16.05190 1.885790
## ACF1
## Training set 0.3126771
## Test set NA
Modeli ETS ka një MAPE të lartë në setin e testimit (16.05%) dhe një RMSE më të madh në krahasim me Holt dhe SES. Kjo mund të sugjerojë se modeli është shumë kompleks dhe nuk është gjithmonë më i përshtatshëm për të gjitha llojet e serive kohore, duke pasur parasysh se mund të ndodhin të dhëna të paqëndrueshme
5.6.5) Vizualizimi i parashikimeve me këto modele
par(mfrow = c(2,2))
plot(model1)
plot(model2)
plot(model3)
plot(model4)
5.7.1) Ndërtimi i modelit
tbats.12<-tbats(ts(data_train$ProdhimiDritherave,start=1961,frequency=1))
tbats.f12<-forecast(tbats.12,12)
autoplot(tbats.f12)+ylab("Prodhimi i Dritherave")+
autolayer(ts(data_test$ProdhimiDritherave,start=c(2009,1),frequency = 1),series="Data")
Modeli TBATS u përdor për të ndërtuar një model të serive kohore mbi prodhimin e drithërave. Modeli i ndërtuar u parashikua për 12 periudha (vite), dhe rezultatet u krahasuan me të dhënat reale të testit.Shikojmë një përshtatje të mirë të vijave.
5.7.2) Rezultatet dhe saktësia e modelit
summary(tbats.12)
## Length Class Mode
## lambda 1 -none- numeric
## alpha 1 -none- numeric
## beta 1 -none- numeric
## damping.parameter 1 -none- numeric
## gamma.values 0 -none- NULL
## ar.coefficients 0 -none- NULL
## ma.coefficients 0 -none- NULL
## likelihood 1 -none- numeric
## optim.return.code 1 -none- numeric
## variance 1 -none- numeric
## AIC 1 -none- numeric
## parameters 2 -none- list
## seed.states 2 -none- numeric
## fitted.values 48 ts numeric
## errors 48 ts numeric
## x 96 -none- numeric
## seasonal.periods 0 -none- NULL
## y 48 ts numeric
## call 2 -none- call
## series 1 -none- character
## method 1 -none- character
accuracy(forecast(tbats.12,12),data_test$ProdhimiDritherave)
## ME RMSE MAE MPE MAPE MASE ACF1
## Training set 5030.823 91887.58 60236.44 -1.300912 10.365632 1.024468 0.1295127
## Test set 21770.164 33000.51 28076.41 3.106058 4.044544 0.477508 NA
Përdorueshmëria e modelit për saktësi: MAPE = 4.04%, që tregon një gabim relativisht të vogël dhe model i përshtatshëm për këtë dataset.
5.7.3) Kontrolli i mbetjeve
checkresiduals(tbats.12$errors)
##
## Ljung-Box test
##
## data: Residuals
## Q* = 7.1862, df = 10, p-value = 0.7078
##
## Model df: 0. Total lags used: 10
Testi Ljung-Box: p-value = 0.7078, që tregon mungesë të korrelacionit të mbetjeve. Testi Shapiro-Wilk: p-value = 6.908e-05, duke treguar se mbetjet nuk ndjekin një shpërndarje perfekte normale. Veme re se autokorrelacionet e mbetjeve jane brenda intervalit te besimit dhe kjo do te thote se seria e mbetjeve eshte stacionare.
5.7.4) Box Test për mbetjet
Box.test(x = tbats.12$errors)
##
## Box-Pierce test
##
## data: tbats.12$errors
## X-squared = 0.16569, df = 1, p-value = 0.684
Një p-value më të madhe se 0.05 (në këtë rast, p = 0.684) që do të thotë se mbetjet janë të pavarura.
5.7.5) Vizualizimi i mbetjeve
library(ggpubr)
# Density plot
ggdensity(tbats.12$errors, fill = "lightgray")
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
# QQplot
ggqqplot(tbats.12$errors)
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
## Don't know how to automatically pick scale for object of type <ts>. Defaulting
## to continuous.
Vizualizimi përfshiu grafikun e densitetit dhe QQ-Plot, duke konfirmuar mungesën e shpërndarjes perfekte normale.Që modeli te jetë i mirë ne presim që kuantilet teorike dhe të zgjedhjes të jenë brenda zonës gri pra afer drejtëzës. Në këtë rast shumica e tyre janë brenda zonës gri.
5.7.6) Shaqpiro Test për mbetjet
shapiro.test(tbats.12$errors)
##
## Shapiro-Wilk normality test
##
## data: tbats.12$errors
## W = 0.86829, p-value = 6.908e-05
Një p-value shumë të vogël (në këtë rast, p = 6.908e-05) sugjeron që nuk mund të pranohet hipoteza e shpërndarjes normale.
Tendencat dhe Analiza e Serive Kohore: • Prodhimi i drithërave ka ndjekur një trend rritës në përgjithësi gjatë periudhës 1961-2020. • Variablat si përdorimi i plehrave kimike dhe emetimet e CO2 treguan ndikim të rëndësishëm në rritjen e prodhimit.
Korrelacioni midis Variablave: • Matrica e korrelacionit tregoi një lidhje të fortë pozitive midis plehrave kimike dhe prodhimit të drithërave (r = 0.89), si dhe midis CO2 dhe prodhimit të drithërave (r = 0.75). • Plehrat kimike rezultuan faktori kryesor që ndikon në prodhimin e drithërave.
Stacionariteti dhe Përgatitja e Modeleve: • Seria kohore për prodhimin e drithërave ishte jo-stacionare dhe u transformua përmes diferencimit për t’u bërë e përshtatshme për modelim. • Modelet u ndërtuan dhe krahasuan për të identifikuar qasjen më të përshtatshme për parashikime.
Modelet e Zbatuara: • Modeli i Regresit Linear tregoi ndikim të qartë të plehrave kimike në prodhimin e drithërave, me një R² = 0.79. • Modelet ARIMA, ETS, Holt-Winters dhe TBATS u përdorën për parashikime afatshkurtra dhe afatgjata. • Modeli TBATS tregoi performancë të mirë për shkak të fleksibilitetit të tij në trajtimin e trendit dhe sezonalitetit.
Performanca e Modeleve: • TBATS rezultoi më i saktë për dataset-in aktual me një gabim relativisht të ulët (MAPE = 4.04%). • Modelet si ETS dhe Holt shfaqën gabime më të larta krahasuar me TBATS dhe ARIMA.
Rezultatet Kryesore të Parashikimit: • Parashikimet treguan se prodhimi i drithërave do të vazhdojë të rritet, por rritja është e ndikuar ndjeshëm nga përdorimi i plehrave kimike. • Në rast të uljes së përdorimit të plehrave kimike ose emëtimeve të CO2, prodhimi mund të ketë ulje të konsiderueshme.
Vlerësimi i Gabimeve dhe Mbetjeve: • Analiza e mbetjeve tregoi që shumica e modeleve nuk kanë autokorrelacion të mbetjeve, duke treguar përshtatje të mirë me të dhënat. • Megjithatë, disa shpërndarje të mbetjeve devijuan nga normaliteti, duke sugjeruar nevojën për përmirësime në trajtimin e të dhënave ose përzgjedhjen e modeleve.
Detyra tregoi se modelet TBATS dhe ARIMA janë të përshtatshme për parashikimin e prodhimit të drithërave në bazë të serive kohore. Përmirësimi i vazhdueshëm i metodave dhe analizave mund të ndihmojë në marrjen e parashikimeve më të sakta për të ndihmuar politikat bujqësore dhe menaxhimin mjedisor.